Attractively bound pairs of atoms in the Bose-Hubbard model and antiferromagnetism 
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We consider a periodic lattice loaded with pairs of bosonic atoms tightly bound to each other via 
strong attractive on-site interaction that exceeds the inter-site tunneling rate. An ensemble of such 
lattice-dimers is accurately described by an effective Hamiltonian of hard core bosons with strong 
nearest-neighbor repulsion which is equivalent to the XXZ model with Ising-like anisotropy. We 
calculate the ground-state phase diagram for a one-dimensional system which exhibits incompressible 
phases, corresponding to an empty and a fully filled lattice (ferromagnetic phases) and a half-filled 
alternating density crystal (anti-ferromagnetic phase), separated from each other by compressible 
phases. In a finite lattice the compressible phases show characteristic oscillatory modulations on 
top of the anti-ferromagnetic density profile and in density-density correlations. We derive a kink 
model which provides simple quantitative explanation of these features. To describe the long- 
range correlations of the system we employ the Luttinger liquid theory with the relevant Luttinger 
parameter K obtained exactly using the Bethe Ansatz solution. We calculate the density-density 
as well as first-order correlations and find excellent agreement with numerical results obtained with 
density matrix renormalization group (DMRG) methods. We also present a perturbative treatment 
of the system in higher dimensions. 

PACS numbers: 37.10.Jk, 03.75.Lm, 05.30.Jp, 75.10.Jm 



I. INTRODUCTION 

Various idealized models describing many-body quan- 
tum systems on a lattice, such as the Heisenberg spin and 
Hubbard models, have been widely studied for decades 
in condensed matter physics [1, 2]. With the recent 
progress in cooling and trapping bosonic and fermionic 
atoms in optical lattices [3] , some of these models can now 
be realized in laboratory with unprecedented accuracy — 
the Hubbard model being a case in point [4]. Imple- 
menting more general models, e.g., extended Hubbard 
or asymmetric spin models, with atoms in optical lat- 
tice potentials is, however, more challenging but poten- 
tially very rewarding. The purpose of the present paper 
is to study an experimentally relevant situation realizing 
the extended Hubbard model or, equivalently, an anti- 
ferromagnetic XXZ model in the Ising-like phase with 
cold neutral atoms in a deep optical lattice potential. 

We consider an optical lattice realization of the Bose- 
Hubbard model with strong on-site attractive interaction 
between the atoms. Specifically, we study a situation 
when each site of the lattice is loaded with cither zero 
or two atoms. Experimentally, this can be accomplished 
by adiabatically dissociating a pure sample of Feshbach 
molecules in a lattice with at most one molecule per site 
[5, 6]. The on-site attractive interaction then results in 
the formation of attractively-bound atom pairs [7, 8] — 
"dimers" , — whose repulsive analog was realized in a re- 
cent experiment [6]. 

For strong atom-atom interaction, either attraction or 
repulsion, the dimer constituents are well co-localized 
[8], and an ensemble of such dimers in a lattice can be 



accurately described by an effective Hamiltonian which 
has the form of a spin-i XXZ model with Ising-like 
anisotropy. The derivation of the effective Hamiltonian 
is given in [9], where we have also discussed its prop- 
erties for the case of repulsive atom-atom interactions. 
Since the resulting nearest-neighbor attraction of dimers 
dominates the kinetic energy, it causes the formation of 
minimal surface "droplets" of dimers on a lattice below a 
critical temperature. In the case of attractive atom-atom 
interaction considered here, the interaction between the 
nearest neighbor dimers is a strong repulsion. We then 
find that the ground state of the system of dimers in a 
grand canonical ensemble exhibits incompressible phases, 
corresponding to an empty and a fully filled lattice as well 
as a half-filled alternating density crystal. These phases 
are separated from each other by compressible phases. 

We calculate numerically and analytically the ground 
state phase diagram for this system in one dimension 
(ID). The critical points can be obtained with the help 
of the Bethe Ansatz making use of the correspondence 
to the XXZ model [10]. In a finite lattice and close to 
half filling, the compressible phases show characteristic 
oscillatory modulations on top of the anti-ferromagnetic 
density profile. A simple kink model is derived which 
explains the density profiles as well as number-number 
correlations in the compressible phases. The long-range 
correlations of the dimer system show a Luttinger liquid 
behavior. We calculate the amplitude and density corre- 
lations in a finite system from a field theoretical model, 
which show excellent agreement with the numerical data. 
The corresponding Luttinger parameter is obtained by 
solving the Bethe integral equations. We then proceed 
to a perturbative treatment of the system in higher di- 
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mensions. Finally, we briefly discuss the implications of 
tunable nearest- neighbour interactions. 

II. EFFECTIVE DIMER MODEL 

We consider attractively-bound dimers on a d- 
dimensional isotropic lattice. Because of the strong on- 
site atom-atom interaction [/ < 0, it is energetically 
impossible to break the dimers, which effectively play 
the role of hard core bosons on the lattice. Via a sec- 
ond order process in the original atom hopping J, the 
dimers can tunnel to neighboring sites with the rate 
J = ~2,P /U > and carry nearest neighbor interaction 
fixed at 4 J. The effective Hamiltonian for the system has 
been derived in [9], 

Hcs = ^ {2ej+U-2dj)mj-j'^cjci+4j'^mjrhi, 

3 (oA {hi) 

(1) 

where and Cj are the creation and annihilation opera- 
tors and rhj = CjCi is the number operator for a dimer at 
site j. In the first term of Eq. (1), the local potential en- 
ergy 2ej of the pair of atoms is modified by an additional 
"internal energy" of the dimer (U — 2dJ^ , which is nega- 
tive for attractive interactions, so that the effective local 
chemical potential is given hy fij = \U\ + 2d J — 2ej. The 
kinetic energy of one dimer described by the second term 
of Eq. (1) spans the interval [—2d J, 2d J] corresponding to 
a Bloch band of a d dimensional square lattice. In com- 
parison, bringing a pair of dimers to neighboring sites 
requires an energy of 8J due to the strong repulsive in- 
teraction in the last term. 

Since the dimers are effectively hard-core bosons, it 
is possible to map the above Hamiltonian onto a spin 
system. Mapping between bosons and spin operators is 
given by the well known Holstein Primakoff transforma- 
tion 

= rhj - 1/2, 
S+ - c] v/r^, (2) 

which preserves the SU(2) commutation relations exactly. 
Since double occupancy is forbidden, rrij — or 1, the 
factor -y/ 1 — rhj is zero for an occupied site and unity for 
an empty site. Therefore, we simply have = and 
Sj' — Cj, so that the equivalent Hamiltonian is given by 

^spi„/J = -J2 {^^^i - + E ^^S-, (3) 

with an effective field of hj — {2ej + U)/ J + 6d. This 
is the anti-ferromagnetic XX Z spin model with a fixed 
anisotropy of 4, i.e., the model is in the gapped Ising-like 
phase. A given total number N of dimers in a lattice of 
Ng sites corresponds, in the spin model, to a fixed total 
magnetization = N — ^Ng. 
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FIG. 1: (Color online) Density of dimers in a ID lattice 
and additional harmonic confinement potential obtained from 
DMRG simulation, with Hj = 18.5 J-2ej and ej / J = jV4400. 
One clearly identifies the incompressible phase with homoge- 
neous filling of (rhj) = 1 in the trap center, and two AF 
phases, separated by compressible intermediate regions. 

III. DIMER SYSTEM IN ONE DIMENSION 

It is now clear that the behaviour of the dimer system 
in one dimension can be determined via isomorphic map- 
ping of Eq. (1) onto the ID integrable XX Z model in a 
uniform field, 

(4) 

where A is the anisotropy of the spin-spin interaction. 



A. Ground-state phase diagram 

An important general feature of the dimer model in 
Eq. (1) is that the ratio of interaction to kinetic energy 
has a fixed value larger than one. As a consequence, 
the ground-state of the system is interaction dominated 
giving rise to interesting correlation properties. 

In a homogeneous system, the ground-state of the sys- 
tem depends only the overall chemical potential /^/J. 
The corresponding phase diagram can be completely 
mapped out in an experiment by adding a shallow ex- 
ternal trapping potential with sufhciently small confine- 
ment such that the local density approximation is valid 
and fj, = \U\ + 2J — 2ej does not change significantly over 
many lattice sites. Then different regions in the trap 
would correspond to different chemical potentials. 

In Fig. 1 we plot the density of dimers in a one- 
dimensional lattice and an additional harmonic trapping 
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potential obtained by numerical DMRG calculations [11]. 
One clearly recognizes three types of regions: In the trap 
center, where the local chemical potential is largest, there 
is a unit filling of dimers. Separated by a spatial region of 
monotonously decreasing average filling follows a region 
where the filling is exactly one half and the dimers form 
a periodic pattern with period 2 and almost maximum 
modulation depth. In this region, the dominant effect is 
the nearest neighbour repulsion 4,/ > of Eq. (1). To- 
wards the edge of the dimer cloud, the average density 
decreases again monotonously to zero. In terms of the 
equivalent spin system, the central region corresponds to 
a gapped phase of full spin polarization caused by a large 
negative effective magnetic field. The region of exactly 
one half average filling corresponds to another gapped 
phase with anti-ferromagnetic (AF) order induced by the 
strong, Ising-like interaction —4:SjSj^i of Eq. (3). The 
intermediate regions are compressible. 

The critical values of the chemical potential for 
the transitions between compressible and incompressible 
phases in ID are known from the work of Yang and Yang 
[10] on the XX Z model of Eq. (4). For the parameters 
of the present system, we have 



IIAF-/J 
fJ-AF+fJ 



-2, (5a) 

8-2715 y ^ 

n^oo ^'^^^ l*^ arccosh(4) j 
3.68361 . . . , (5b) 

(-1)" 



8+2VT5 



cosh (n arccosh(4)) 



n=—oo 

12.31638..., (5c) 
18. (5d) 



These values agree very well with those obtained from ex- 
act diagonalization on a small homogeneous lattice with 
Ng = 10 sites and periodic boimdary conditions, as well 
as DMRG simulation with up to Ns = 300 and open 
boundary conditions. They also match the different re- 
gions of Fig. 1. 



B. Mott-insulating phcises 

In the language of spin Hamiltonian, phases with zero 
(A^ = 0) or full {N = Ns) filling c;orrespond to ferromag- 
netic phases with a simple form of the ground state 



IV't) = |T,T,T,...,T). 



(6) 
(7) 



Particle-hole excitations are not possible in the Mott- 
insulating state (7), while inserting a particle into (6) or 
removing one from (7), corresponding to flipping a spin, 
carries finite energy cost given by Eqs. (5a) and (5d). 
Hence these phases are incompressible. 



For half filling (N = ^Ng) the situation corresponds 
most closely to an AF phase. However, in this case the 
simple Neel state 



l<^)= |...,i,T,i,T,i,T,i,T,...), 



(8) 



is not an exact eigenstate of the full Hamiltonian Hcs 
m (1). Rather, [V^i"^) is an eigenstate of Hamiltonian 



H. 



(0) 
eft 



HeS — i^hop without the hopping term Hh, 



op 



— J^^(c]_,_^Cj+cjcj+i). Due to Hhop a dimer can tunnel 
from an occupied site to a neighboring empty site, which 
in terms of the Neel state (8), corresponds to flipping two 
neighboring spins, resulting in a state of the form 



Tj-i-i) T) i> T, • • •)• 



(9) 



If we assume periodic boundary conditions and an even 
number of lattice sites Ng, there are j = 1, . . . , Ng differ- 
ent states (9), one for each link where two neighboring 
spins can be flipped. Each of those states IV'j^'') has a 
larger repulsive (Ising) interaction energy which is 

increased by 8 J relative to energy E^p of state IV'af)- 
is tempting to treat the the smaller hopping -ffhop as per- 
turbation with respect to H^g , but unfortunately already 
the first order correction carries a contribution from all 
Ng possible states in Eq. (9). In higher order perturba- 
tion theory the number of contributing states increases 
with higher powers in Ng, so that the perturbation series 
diverges in the thermodynamic limit. 

However, if we are interested in local observables, such 
as the density in Fig. 1, it is possible to restrict the per- 
turbation only to those hopping terms which change the 
value of the density at a particular point. In particular, 
in order to calculate the ground state expectation value 
of (■0af|5'||^af), we have to make the following ansatz 
for the ground state |i/'af) 



E' 



(0) 
AF 



E. 



(1) 



(10) 



which can be normalized by a factor of l/^/l -I- 1/32. 
This state is in general a bad approximation to the 
ground state, but it describes very well which terms in 
the Hamiltonian affect the local density at site j, since 
higher order hopping only contributes 1/64 or less. Ac- 
cordingly, the local density is given by 



(VafI^II^af) 



^ ^ 66 V 32 



(-ly— (11) 



which corresponds to a deviation of about 0.03 from per- 
fect alternating order. Our numerical results for homoge- 
neous systems show a deviation of about 0.032 which is in 
very good agreement with the prediction. Even though 
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the half filling state always implicitly contains excitations 
of type (9), the removal or addition of a particle still costs 
relatively large energy given by Eq. (5b) or (5c), which 
makes the AF phase incompressible. 

C. Properties of compressible phases 

In the remainder of this Section, we examine the com- 
pressible phases, mainly in the vicinity of the antiferro- 
magnetic phase, using two different approaches. The first 
is perturbative in nature and relies on the fact that the 
nearest neighbor interaction energy between the dimers 
exceed the dimer hopping energy by a large factor of 8. 
We show that the system can approximately be treated as 
a non-interacting gas of kinks that behave like hard-core 
bosons. The second approach aims to describe long-range 
correlations employing the Luttinger-liquid theory. The 
relevant Luttinger parameter can be obtained by Bethe 
Ansatz solution of the equivalent XX Z spin model. 

1. Non-interacting kink approximation 

In Fig. 2 we plot the density distribution of dimers in 
a homogeneous lattice of Ng = 99 sites obtained from 
DMRG simulations with different number of dimers N. 
An infinite (hard-wall) confining potential eo — eioo 
+O0 has been used, which imposes on Hamiltonian (1) 
open boundary conditions with toq — ^tiioo = 0. Due 
to the asymmetric coupling at the boundaries, the end 
sites j — 1 and j = 99 prefer to be occupied with a 
particle. To accommodate an oscillating density wave, we 
therefore use odd number of lattice sites Ns- Note that 
the open boundary condition for the particles in Eq. (1) 
corresponds to an additional effective edge field hi — 
hgg = —2 for the spins in the XXZ-mode\ of Eq. (3), 
which has the analogous effect of polarizing both end 
spins up. 

As seen in Fig. 2, the ground state for = 50 ex- 
hibits density oscillations corresponding to the AF Neel 
order, up to the small correction discussed in Sec. IIIB. 
Adding particles leads to modulated density distribution, 
with the envelope of modulation having regularly spaced 
nodes whose number is equal to twice the number of addi- 
tional particles. In the following we will provide a simple 
theoretical understanding for this effect. 

Without the small hopping term i/hop, the ground 
state of Hamiltonian (1) for half filling is the AF state 

IV'af) of (^)' which is twofold degenerate. The AF 
order with period 2 effectively doubles the size of the 
unit cell. Adding then a particle to IV' af) costs exactly 
an energy of {h + 8) J, resulting in state 

I ■ • ■ J T: Tj: Ti + ll Ti Tl il T: ■ ■ 

which is energy-degenerate with any state of the form 

= |...,T,i,T„T,i,...,i,T,',T,i,...>. (12) 
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FIG. 2: (Color online) Particle density profile in a homoge- 
neous lattice with Ns = 99 sites and hard-wall boundaries, 
for different particle number A''. For half filling, A'^ = 50, 
the ground state has nearly perfect AF order. Adding one, 
two and three particles leads to the density-wave modula- 
tions with the number of nodes equal twice the number of 
additional particles. 



Hence, the additional particle causes effective domain 
walls, which can be placed anywhere in the system and 
play the role of mobile kinks at positions j and j' between 
AF regions with different orientation. Note that without 
hopping any number of particles above half filling can 
be created at the critical field hi ■* = — 8 and placed in 
an arbitrary arrangement as long as no two neighbor- 
ing lattice sites are empty. In other words, at hi ^ we 
have a huge degeneracy of states with any magnetiza- 
tion AI2 > corresponding to arbitrary arrangement of 
antiferromagnetic regions and spin-up ferromagnetic re- 
gions. The analogous statement is also true at the upper 
critical field h[^^ = 8, where the degenerate subspace is 
defined as states with < where no two neighboring 
spins may point up. This degeneracy implies that with- 
out hopping the transition from the antiferromagnetic in- 
compressible phase to the ferromagnetic incompressible 
phases is infinitely sharp at the effective critical magnetic 
fields hi'^K As we will see below, however, the hopping 
lifts this degeneracy and therefore is crucial for the sta- 
bility of compressible phases over finite ranges of field h 
as observed in Fig. 1. 

The hopping H^op is also responsible for the modu- 
lated wave patterns seen in Fig. 2. Starting from the 
AF state in Eq. (8), we now insert more and more par- 
ticles each producing a pair of kinks. The states in 
Eq. (12) can be considered as AF states with a pair 
of kinks, one at even sites and one at odd sites; e.g., 
state I ii,T2,i3,T4,T5,T6,i7,T8,i9,---) has kinks at sites 
4 and 5, while state | ii, T2, is, T4, T5, ie, Tr, Ts, ig, ■ • ■) 
has kinks at sites 4 and 7. i?hop has non-vanishing ma- 
trix elements within the subspace of energy degenerate 
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FIG. 3: (Color online) Top: ID chain with one particle added 
to the AF state creating a pair of odd (red) and even (blue) 
kinks. The hopping Hamiltonian Hhop leads to a motion of 
the odd and even kinks on odd or even sites, respectively. In- 
terchange of odd- and even-site kinks is not possible. Bottom: 
mapping onto an effective lattice with lattice constant 2. 



states with fixed number of additional particles. Within 
this manifold of states, hopping of the additional particle 
corresponds to free motion of kinks, wherein an even-site 
kink moves only on even sites and an odd-site kink on odd 
sites, as illustrated in the top part of Fig. 3. Furthermore, 
the even and odd site chain kinks cannot exchange their 
relative order. Note that hoping of a particle surrounded 
by two empty sites is energetically suppressed. 

Given a fixed number of additional particles or holes, 
the motion of the corresponding kinks is equivalent to 
the motion of hard-core bosons in an effective lattice with 
lattice constant 2. To sec this, consider the case of q addi- 
tional particles on top of the half filled lattice; the oppo- 
site case of holes follows from the particle-hole symmetry. 
Let the positions of the kinks be ji < j2 < ■ ■ ■ < j2q- If 
ji is even (odd) then is, js, jV, . . . are also even (odd) and 
j2, j4,i6: ■ • ■ are odd (even). We now perform a mapping 
onto a new lattice which we call the kink lattice. The 
quasi-position kn of the nth kink is then 



2 

jn+n 



if ji is even 
if ji is odd 



(13) 



This mapping is illustrated in the lower part of Fig. 3. 

Evaluating the matrix elements of the hopping Hamil- 
tonian Hiiop in the subspace of states with constant num- 
ber of kinks, we find that the latter can be treated as 
hard-core bosons or non-interacting fermions on the kink 
lattice, only if we consider the absolute value of the wave- 
function. The corresponding hopping strength on the 
period-2 lattice is again J. The exchange symmetry can- 
not be determined straightforwardly and therefore we 
employ this approximation only to determine the den- 
sity distribution of dimcrs. For simplicity wc choose the 
fermionic exchange symmetry. 

Let us assume that the lattice is large and consider a 
particle filling close to the antiferromagnetic case. In this 
limit, the kinks can be regarded as moving on a contin- 
uum. This means that the dynamics of the kinks can 
now be determined by solving the Schrodinger equation 
for non-interacting fermions. For N = + 1, i.e.. 



one additional particle, we have a pair of kinks whose 
ground-state wavefunction is 



. TTXi . 2nx2 . 7rX2 . 27ra;i 
sm — — sm — sm — — sm — - — 



W2(xi,a;2) = — 

(14) 

where L = ^Ng -f 1 is the length of the kink lattice. The 
left-most kink shall move on the odd sites. A particle is 
sitting on an even site j if and only if one chain kink is 
to the left of j. Thus the density of particles on the even 
sites is 

{m{x))^2f dyi f dy2*;(yr,y2)*2(yr,2/2). (15) 

The prefactor of two emerges here because the integral 
occurs twice with interchanging the roles of yi and y2- 
Although straightforward, we do not give the analytic 
expression of Eq. (15) since it is rather long. At the odd 
sites we get accordingly 

l~{rh{x)) = dyi dy2-i>2{yi,y2)-^2iyi,y2) 
Jo Jo 

dyi / rfy2 *2(2/i.y2)*2(yr,y2)- (16) 



With q additional particles, the fermionic ground state 
wavefunction for 2q kinks is 



'^2q{xi 



,a;2q) = ^ ^^^= n 0P(„)(a;„), (17) 



where the sum is over all permutations P of numbers 
{1,2,3, ...,2(7} and 



4>n(x) 



2 irnx 



with L = ^Ng+q. This results in the density distribution 



sgn(P) sgn(g) 
{2k + l)\{2q - 2k ~ 1) 



fe=0 P,Q 
2k+l 2q 

X l[ l{0,x,P{n),Q{n)) l{x, L, P{n),Q{n)) 

n=l n=2k+2 

(18) 

where Q denotes the permutations of {1, 2, 3, . . . , 2(7}, 
and 



/(a, 6, 71, m) 



dx(t>l{x)(t>m{x), 



with n,m G {1, 2, 3, ... , 2q}. In Eq. (18) we have taken 
into account that there are (2g-2A;-i)'!'(2/c-n)! possibilities 
of choosing 2fc + 1 kinks to the left of j. 

The dashed red lines in Fig. 2 show the analytic results 
for the particle density in a box potential with the lattice 
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FIG. 4: (Color online) Particle density- density correlations 
in a homogeneous lattice with Ns — 99 sites and hard-wall 
boundaries, for different particle number A''. The blue lines 
correspond to numerical DMRG results, the red dashed lines 
to the predictions of the kink approximation. 



2. Field theoretical approach 

The spin chain equivalent to the dimcr Hamiltonian (1) 
in ID is given by Eq. (4) with A = 4 and open boundary 
conditions. At zero magnetization = 0, the XX Z 
model (4) is gapped, since A > 1. However, as described 
in Sec. Ill A, the gap can be closed by a field between 
the two critical values, hi ^ < h < hi'^\ In other words, 
the system is critical for any nonzero magnetization away 
from the fully magnetized case. In this regime, the lead- 
ing low-energy effective theory is a Luttinger liquid with 
two parameters, the spin velocity v and Luttinger pa- 
rameter K. These are functions of the magnetization 
per site = M^/Ns and anisotropy A [12] (which for 
the particular dimcr model here is fixed, A = 4). 

In order to calculate correlation functions, we first de- 
rive the Luttinger parameter K{s^) from the exact solu- 
tion [13]. We write 



K = e{b), 



(20) 



where the function ^{x) is determined by the integral 
equation 



filling slightly above one half obtained from the kink ap- 
proximation. The agreement with the numerical DMRG 
data is rather good. The kink model also explains in a 
very intuitive way the pairwise appearance of nodes with 
adding every particle to the lattice. 

Particle number correlations can be derived in the 
same manner. For two even sites at positions ji and j2, 
the configurations contributing to the correlations cor- 
respond to an odd number of particles to the left of ji , 
an even number of particles between ji and j2, and an 
even number of particles to the right of j2- The particle 
density-density correlations are then given by 



<(,-!) 



sgn(P) sgn(Q) 



L(2fci + l)!(2fc2)!(2A:3 + l)! 



(m(x)m(y)) = 

fel,fc2,fc3=0 P.Q 

2ki + l 2fei+2fe2 + l 

X n l{0,x,P{n),Q{n)) [] l{x,y, P{n),Q{n)) 



n=l 



n=2fei+2 



2q 



X n l{y,L,P{n),Q{n)) 

n=2fei+2fe2+2 



for X < y. (19) 



In Fig. 4 we plot the density-density correlations ob- 
tained from DMRG calculations (blue solid line) and 
the kink model (dashed red lines), displaying very good 
agreement. We finally note that within the approxima- 
tion of non interacting kinks, first order correlations exist 
only between neighboring sites. This perturbative model 
therefore can not accurately describe such correlations. 



ax) = l+ I n{x-y)ay)dy, (21) 

-b 



with the kernel 



k{x) 



sinh 2r/ 



1 



TT cos 2x — cosh 2?7 



A = coshry > 1. 



The parameter b in Eqs. (20) and (21) is implicitly de- 
fined through 



1 



p{x) dx, 



(22a) 



-6 



p{x) = d{x)+ / K{x~y)p{y)dy, (22b) 



Filling 



FIG. 5: Dependence of the Luttinger parameter K on the 
mean lattice filling N/Ns, for A = 4. 
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where 

,/ ^ 1 sinhr? ^ 
d{x) = — — , A = coshr?>l. 

TT COS 2x — cosh r] 

Equations (20)-(22) are solved numerically by discretiz- 
ing the integral and inverting the resulting matrix equa- 



tion. Figure 5 shows the function K{s^) for A = 4. 

Within the Luttinger liquid approach, one- and two- 
point correlation functions can be calculated using the 
standard mode expansion of bosonic fields [14] for open 
boundary conditions [15, 17]. Then the spin-spin corre- 
lation function in the ground state reads 



{S^{x)S^{y)) = [s^Y-B 



K 



1 



8(iV, + If sin2 JL^ 



-vL 

2(JVs + l) 



2 Tr{x+y) 



cos [{2kF + e/Ns)x + ipi] 



sm' 



2(JV3+1) 



2K 



K 



+C2 



cos[{2kF + 0/N,)y- 



K 



f2] _^ ^ cos [(2fcf + 0/Ns) x + S 



Ns + l 



Try 



K 



7r(:r+;y) 
2(iV,,+l) 

2(W, + 1) 



(23) 



with the Fermi wavevector kp = 7r(l — 2s^)/2. Here the amplitudes B, Ci,2, D, the shift 9, and the phases </?i,2, S result 
from bosonization of operators on the lattice. We consider them as parameters in Eq. (23) that are fixed numerically 
by fitting to the DMRG data. The exponents, however, arc obtained from the Luttinger liquid parameter K , which is 
given by the Bethe Ansatz. Figure 6 shows the remarkable agreement between the two approaches. Note the shift in 
the wavevectors of the oscillations by a constant that depends on the boundary conditions, the interaction and the 
magnetization. It has also been observed in the context of density oscillations in the open Hubbard model [16, 18]. 
The corresponding result for the first-order correlation function in the ground state is 
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cos [{2kF + e/Ns){x-y) + 6] 



^sir 



K 



sm 



sm 



T(x+y) 
2(jV,, + l) 
Tr{x — y) 



2K 



, ^ cos [{2kF + e/N,)x + ipi] „ cos [{2kF + e/Ns)y + <P2] 
—5= 1- O2 — 



K 



(24) 



Similarly to Eq. (23), the quantities B,C\^2,D,5,(f)i^2,9 
are considered as fitting parameters. The resulting curves 
are shown in Fig. 7. 

IV. PHASE DIAGRAM IN HIGHER 
DIMENSIONS 

We now derive the phase boundaries for the dimer sys- 
tem in two and three dimensions. To that end, we employ 
the strong-coupling approach [19], wherein the hopping 
term -ffhop = i) ^j*^* '-'^ Hamiltonian (1) is treated 

as small perturbation with respect to = Hes — Hhop- 
a. Zero-hopping limit. Without the hopping, the 
grand canonical operator for the dimer system reads 

^is = 4 J ^ rn^mi - M X] ■ (^5) 

In this (formal) limit, the model is isomorphic to the 
Ising model in an external magnetic field and it has two 



critical points 

= 0, (26a) 
mJ^V-^ = 16d (26b) 

For very small values of the chemical potential, /i < 0, 
all spins are polarized in the —z direction, which in the 
dimer language corresponds to a state with zero dimers 
at each lattice site. For sufficiently large values of chem- 
ical potential, fj, > 16dJ, all spins are aligned in the +z 
direction, i.e., we have unit filling of the dimer lattice. Fi- 
nally, for intermediate values of the chemical potential, 
< fi < ^^^\ the ground state is twofold degener- 
ate and has antiferromagnetic order corresponding to a 
"checkerboard-crystal" lattice of dimers. 

b. Boundaries of ferromagnetic phases. When the 
hopping term i/hop is brought into the picture, the two 
critical points extend to two critical regions in which the 
system is compressible. In order to determine the chem- 
ical potentials at which the transitions between the com- 
pressible and incompressible ferromagnetic phases take 
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FIG. 6: (Color online) {S%x)S'^ {50)) correlations obtained 
from the DMRG (dots) and the Luttinger-liquid approxima- 
tion (solid lines). 
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FIG. 7: (Color online) {5"''(a;)5'~(50)) correlations obtained 
from the DMRG (dots) and the Luttinger-liquid approxima- 
tion (solid lines). 



tributes 8 J of repulsive interaction energy and we obtain 

E{N = N,) = {U -2dJ)Ns + 8dJN,, 
E{N^N,-l) ^ {U - 2dJ)iNs - 1) 

+8dJ{N, - 2) - 2dJ. 

The critical chemical potentials, are then determined by 
the energy difference E{1) — E(0), at which we add a 
dimer to the empty lattice, and E{Ns) — E{Ns — 1), at 
which we add a hole to (or remove a dimer from) the 
filled lattice 



Mx/J = -2d, 
/^l/J = I8d. 



(27a) 
(27b) 



It should be noted that the hopping Hamiltonian does 
not modify the corresponding states in the two insulat- 
ing phases, i.e., within the empty and fully-filled lattice 
phases there are no fluctuations of the dimer number, 
which is exactly zero or one per site, respectively. 

c. Boundaries of antiferromagnetic phase. We now 
calculate the lower and upper critical chemical potentials 
Maf=f for the AF phase, up to the second order in dimer 
hopping. At exactly half filling, N = ^Ng, the ground 
state is an almost perfect antiferromagnet with an alter- 
nating density structure in which there are no nearest- 
neighbour repulsive junctions between the dimers. How- 
ever, due to Hhop each dimer undergoes highly nonrcs- 
onant transitions to the neighboring empty sites, whose 
number is 2d (see Fig. 8 illustrating the 2D case), and we 
obtain 



E^^^i^Ns) = {U~2dJ)\Ns 



8J(2d- 1) 



\Ns 2d, 



where the last term describes the second-order energy 
shifts resulting from the virtual transitions of the dimers. 

In the cases of = ^A^s ± 1, the added dimer or dimer 
hole can not freely move for d> 2, since it would require 
two hopping events, as can be seen in Fig. 8. This should 
be contrasted with the ID situation, wherein adding a 
dimer or a hole to the AF phase creates two mobile kinks. 



place, we calculate the particle and hole excitation ener- 
gies for zero and full filling of a finite lattice with even 
number of sites A^^ using Hamiltonian (1). 

In the case of a single dimer in an empty lattice, there 
is no contribution from the interaction energy and we find 
immediately without resorting to perturbative treatment 
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FIG. 8: (Color online) AF state in 2D lattice with an addi- 
tional particle (left) or hole (right). Virtual hopping of parti- 
cles (holes) adjacent to the defect (red) and in the bulk (blue) 
lead to different second-order energy contributions. Hopping 
of the additional particle or hole is not allowed. 



discussed in Sec. IIIC 1. Taking into account the second 
order corrections due to virtual transitions of dimers ad- 
jacent to the extra dimer or dimer-hole (see Fig. 8 left or 
right, respectively), we obtain 



{U -2dJ){\Ns ± 1 
{J? 



8J(2d- 1) 

{J? 
8J{2d-2) 



{^Ns - 2d) 2d 
2d{2d- 1). 



The lower haf- = E^^^ {^N,) - E'^^^ (iiV, - l) and up- 
per i^AF+ = -E^^) (iiV, + 1) - E^'^'> (^Ns) critical chemical 
potentials for the AF phase in d > 2 dimensions are then 
given by 



fJ-AF-/J 



d 



4(2d- l)(2d- 2)- 
HAF+fJ = I6d- 



4(2d- l)(2<i- 2)' 



(28a) 
(28b) 
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FIG. 9: Top: Ground state density profile in a homogeneous 
lattice with Ns = 99 sites and hard-wall boundaries, for N = 
50 particles, but with hopping strength equal to the nearest 
neighbor repulsion, A = 1 (compare to Fig. 2, top panel). 
Bottom: The same, but without the effective edge fields in 
the spin chain model. 



V. THE ROLE OF ANISOTROPY A 

The effective Hamiltonian (1) has a fixed relation of 
the nearest neighbour interaction to hopping, which re- 
sults in a fixed Ising like anisotropy A = 4 in Eq. (3). 
It is now interesting to also consider the more general 
case of tunable anisotropy, which could for example be 
realized with dimers consisting of two different atomic 
species [20]. As the hopping becomes stronger, the per- 
turbative analysis used in the previous sections becomes 
unreliable. It is known that the model of Eq. (4) is critical 
for — 1 < A < I. Therefore, the perturbation treatment 
breaks down exactly at the point where the hopping be- 
comes equal to or larger than the nearest neighbor dimer- 
dimer interaction. However, it is still possible to use the 
field theoretical methods of Sec. Ill C 2 to calculate the 
correlation functions and expectation values. 

In the critical region excitations are gapless in the 
thermodynamic limit, so that there is no incompressible 
phase at half-filling. The crossover between the com- 
pletely filled and completely empty regions in Fig. I is 
therefore continuous as a function of the effective field 
and there is no extended half-filled phase. The strength 
of hopping is therefore crucial for the behavior of the 
system: weak hopping enables the presence of a com- 
pressible phase between the incompressible ferromag- 
netic and antiferromagnetic phases. With increasing the 
hopping strength the incompressible antiferromagnetic 
phase shrinks and completely vanishes when the hop- 
ping reaches the value of the nearest neighbor interaction, 
A < 1. 

Equally interesting is the effect of hopping on the den- 
sity of dimers along the chain. As discussed in Sec. Ill C, 
density- wave modulations appear in Fig. 2 because of the 
effective motion of kinks. The complex interplay between 



the kinetic and interaction terms in the critical region 
now leads to a further modification of the density pattern 
along the chain [21-23]. In particular, the amplitude of 
the ground state density oscillations is now significantly 
reduced towards the middle of the chain with a charac- 
teristic drop-off as shown in Fig. 9. Excited states with 
larger N would then exhibit modulations on top of this 
ground state pattern, similar to Fig. 2. 

Interestingly, the exact form of the boundary condi- 
tions play now a much more important role. Namely, the 
effective edge field of the spin chain model discussed in 
Sec III C 1 accounts for a large part of the ground state 
density oscillations. For comparison, in Fig. 9 we also 
show the density for the spin chain model without any 
edge field. The density amplitude is now smaller near 
the boundary. At A = 1, this amplitude has been pre- 
dicted to follow approximately a y^sm{7Tx/Ns) behavior 
[21], which however is strongly affected by temperature 
[22, 23] due to the gapless modes. 

VI. SUMMARY 

In this paper, we have studied the many-body dynam- 
ics of attractively bound pairs of atoms in the Bose- 
Hubbard model. When the on-site interaction between 
the atoms exceeds by a sufficient amount the bandwidth 
of the lowest single-particle Bloch band, the pairs are well 
co-localized and can be treated as composite dimer par- 
ticles. Then the effective model for dimers on a lattice 
is equivalent to an asymmetric spin-^ XX Z model in an 
external magnetic field: the nearest neighbor interaction 
between the dimers translates into an Ising-type spin- 
spin interaction and the dimer tunneling to a spin-spin 
coupling in the x — y plane. 
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The case of repulsively bound pairs studied in [9] cor- 
responds to a ferromagnetic Ising coupling. In contrast, 
for attractively bound pairs analyzed here, the Ising cou- 
pling is antiferromagnetic leading to a much richer phase 
diagram. The asymmetry parameter A of the XXZ 
model is equal to 4; as a result, the system is gapped for 
both zero and full magnetization. The zero magnetiza- 
tion state, corresponding to exactly half filling of dimers, 
exhibits antiferromagnetic order, while the full magneti- 
zation (ferromagnetic) states correspond to vanishing or 
full filling of dimers. When the effective magnetic field, 
or for that matter the chemical potential for the dimers, 
exceeds critical values, the gap is closed and the dimer 
system becomes critical with finite compressibility. 

We have derived the critical values of the chemical po- 
tential for the transition points from the gapped ferro- 
magnetic and antiferromagnetic phases to the compress- 
ible phases in ID, employing the known exact solutions 
of the XXZ model, and by using a strong-coupling ap- 
proximation in higher dimensions. The properties of the 
compressible phases are quite different in one and higher 
spatial dimensions. In ID, close to half filling, the system 
can be well described by kink-like domain walls which 
separate antiferromagnetic strings of opposite phase and 
can propagate through the lattice almost freely. In higher 
dimensions, the motion of similar defects is strongly sup- 



pressed. A simple approximate description in terms of 
non-interacting kinks gives rather accurate predictions 
for the dimer density as well as non-local density-density 
correlations. 

In order to explain the first order correlations, we em- 
ployed a field theoretical approach based on the Lut- 
tinger liquid theory. The corresponding Luttinger pa- 
rameter was obtained by solving the Bethe Ansatz equa- 
tions for the equivalent XXZ model in the regime of crit- 
ical magnetic fields. The expressions for the first-order 
and density-density correlations showed remarkably good 
agreement with the numerical data obtained by DMRG 
simulations. Finally we discussed the consequences of 
changing the anisotropy parameter of the XXZ model. 
Our studies attest that interaction bound pairs of atoms 
in deep optical lattices can provide a versatile tool to 
simulate and explore quantum spin models. 
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